Collapse and revival oscillations as a probe 
for the tunneling amplitude in an ultra-cold Bose gas 
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We present a theoretical study of the quantum corrections to the revival time due to finite tun- 
neling in the collapse and revival of matter wave interference after a quantum quench. We study 
hard-core bosons in a superlattice potential and the Bose-Hubbard model by means of exact numer- 
ical approaches and mean-field theory. We consider systems without and with a trapping potential 
present. We show that the quantum corrections to the revival time can be used to accurately deter- 
mine the value of the hopping parameter in experiments with ultracold bosons in optical lattices. 
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I. INTRODUCTION 

Collapse and revival oscillations in a Bose-Einstein 
condensate loaded in an optical lattice were first exper- 
imentally observed in 2002 [l[ and since then have been 
a subject of much theoretical and experimental interest 
. This phenomenon is understood as an oscillation be- 
tween an initial coherent state and a final non-coherent 
(collapsed) state in a lattice where, after a quench, the 
hopping parameter between sites is negligible. Very re- 
cently it has been argued that such collapse and revival 
oscillations can be used as a very sensitive probe for ef- 
fective three-body and higher interactions @ by studying 
the time evolution of the visibility of the interference pat 
tern. This has been investigated theoretically in Ref. 
It was assumed both in experiment and in theory [H, 
that the initial state is a coherent state and that after 
the quench the tunneling amplitude is negligible, that 
is, that the systems were in the atomic limit. Following 
these assumptions, the time-evolving state is a product 
of coherent states localized at each lattice site and one 
deals with an effective one-site problem. 

Our goal in this article is to go beyond the previ- 
ous analysis and present a full many-body study of col- 
lapse and revival phenomena in one-, two-, and three- 
dimensional cubic lattices. We study hard-core bosons 
in the presence of a superlattice and the Bose-Hubbard 
model. For both cases we consider systems without and 
with a trapping potential present. For the hard-core case 
we use exact numerical approaches in one and two di- 
mensions and compare them with the predictions of a 
Gutzwiller mean-field theory. We show that the latter is 
qualitatively and quantitatively correct when determin- 
ing the revival time for small hopping amplitudes. Build- 
ing on that, we present an analysis for the Bose-Hubbard 
model that is solely based on the Gutzwiller mean-field 
approach. We also provide an analytical solution for the 
homogeneous hard-core case. 

A previous study Q considered the effect of a finite 
hopping on the damping of the collapse and revival os- 
cillations in a lattice without a confining potential. The 



quantitative results presented there applied to an initial 
coherent state as in the articles mentioned previously. 
In contrast, we study the dynamics starting from initial 
states that are the exact many-body ground state in some 
cases and the appropriate Gutzwiller ansatz in the other 



We find the functional form of those corrections to 
the revival time in the atomic limit, and show that if 
one knows the values of the superlattice potential for the 
hard-core case or the onsite interaction U for the Bose- 
Hubbard model, such corrections can be used to accu- 
rately determine the tunneling amplitudes in experiments 
in optical lattices. In the atomic limit, for the Bose- 
Hubbard model, the revival time t rev is t^° m = 2irh/U 
[H, H| . Our general strategy is to numerically calculate 
the deviations from this value for < J/U < 1. Given 
that in experiments with ultracold gases in optical lat- 
tices the hopping parameter J is exponentially sensitive 
to the lattice depth, while the onsite repulsion U is only 
power-law dependent on the lattice depth, this method 
of determining J by means of the revival time may be 
more accurate than the approaches followed so far. Our 
results are also relevant to cases in which the tunneling 
is small but one is still interested in estimating its value. 

The article is organized as follows. In Sec. UH we study 
the collapse and revival in the hard-core case in the pres- 
ence of a superlattice potential. This is done using nu- 
merically exact methods in one and two dimensions. In 
Sec. IIII1 we introduce the time-dependent mean-field ap- 
proach to the hard-core boson problem and compare its 
results with numerically exact ones in order to quantify 
the predictive power of the mean-field approximation for 
the revival time. In addition, we make some general theo- 
retical statements and provide an analytical solution for 
the case without the trapping potential. Section IIVI is 
devoted to analyzing the Bose-Hubbard model in three 
dimensional systems of soft-core bosons with and with- 
out confining potentials present. Finally, we present our 
conclusions in Sec. IVl 
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II. NUMERICALLY EXACT RESULTS FOR 
HARD-CORE BOSONS IN A SUPERLATTICE 

We first introduce the model Hamiltonian considered 
to study hard-core bosons in a superlattice. We also dis- 
cuss its relation to the interaction quench in the soft-core 
case that is studied in Sec. IIVI 



A. Motivation and methods 

The hard-core boson Hamiltonian on a superlattice 
with period 2 in the presence of a harmonic confining 
potential can be written as 
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(1) 

where the hard-core creation and annihilation operators 
at site i are denoted by b\ and b i , respectively, and the lo- 
cal density operator by hi = Sj^. Hard-core boson oper- 
ators satisfy standard commutation relations for bosons. 
However, on the same site, they also satisfy the constraint 
b\ 2 =bf = 0, which precludes multiple occupancy of the 
lattice sites. The other parameters in Eq. ([1} are the 
hopping constant J between nearest-neighbor sites {ij), 
the strength of the superlattice potential A, and the cur- 
vature of the harmonic trap V. The distance from site 
i to the center of the trap is measured in units of the 
lattice constant a, which we set to unity. In what follows, 
we also denote the total number of lattice sites by L and 
the total number of particles by N. 

The hard-core model on a superlattice is particularly 
suitable to study collapse and revival phenomena be- 
cause in one dimension (ID) it can be exactly solved 
by means of the Jordan- Wigner mapping to nonintcr- 
acting fermions [16[. In equilibrium, these systems were 
studied in Ref. |8[, where they were shown to exhibit 
ground state phases that are similar to those of the Bose- 
Hubbard model. The nonequilibrium dynamics of hard- 
core bosons in a superlattice potential was studied in Ref. 
Q , where collapse and revival oscillations of the momen- 
tum distribution function were observed. 

A quench of the superlattice potential A in Eq. ([T]) has 
a similar effect to a quench of U in the Bose-Hubbard 
model. From a simple band-structure calculation it fol- 
lows that in ID A opens a gap T = 2A in the hard-core 
boson energy spectrum [H, Q , 



e±(fc) = ±sJU 2 cos 2 {ka) + A 2 , 



(2) 



where "+" ("— ") denotes the upper (lower) band. In two 
(2D) and three (3D) dimensions, hard-core bosons cannot 
be mapped to noninteracting fermions. The phase dia- 
grams for the ground state of such systems were studied 
in detail in Refs. (Tol . ITl| by various numerical and an- 
alytical approaches, and were shown to be qualitatively 



similar to the phase diagrams of the Bose-Hubbard model 
in 2D and 3D PJ-G3. 

As already noted, in the atomic limit of the Bose- 
Hubbard model, the revival time after the interaction 
quench is given by t^° m = 2ir/U (we set ft = 1 hence- 
forth); similarly, for hard-core bosons in a superlattice 
potential, it follows that t^° m = tt/A [|. 

Hence, in this section we take advantage of the fact 
that the nonequilibrium dynamics of hard-core bosons in 
ID can be exactly solved for large system sizes by means 
of the approach presented in Refs. [15| . which makes use 
of the Jordan- Wigner transformation to noninteracting 
fermions [Hj]. In 2D, due to the reduced Hilbert space 
(when compared to soft-core bosons), one can perform 
full diagonalization calculations for small, but meaning- 
ful, periodic systems. All our exact results for 2D hard- 
core systems were obtained in 4 x 4 lattices with periodic 
boundary conditions. 

The two preceding approaches allow us to make ex- 
act predictions for the quantum corrections due to finite 
hopping amplitudes to the revival time ij! 
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which in turn will help us gauge the accuracy of the mean- 
field approach that we use later for studying the Bose- 
Hubbard model. 

In the latest experimental and theoretical studies 
the main observable under consideration was the visibil- 
ity of the interference pattern. Here, instead, we focus 
our attention on the time evolution of the n&=o momen- 
tum occupation number, 



n k =o 



(4) 



which is also measured in experiments. 



B. Results for hard-core bosons in a periodic 
potential in ID and 2D 

In all cases for hard-core bosons presented here we con- 
sider the following quench. A system is prepared in a su- 
pcrfluid state with J lnl = 1 (which sets our energy scale) 
and no superlattice potential. At time t = the superlat- 
tice potential A is quenched to a constant value Afi n = 1 
and the hopping constant J is reduced to several val- 
ues Jfi n < 1 (in the remainder of the text, the notation 
Afi n = A and Jfi n = J is used in all unambiguous cases). 

In Fig. [TJ we show the time evolution of nk=o after 
this quench in (a) a chain and (b) a 4 x 4 cluster, both 
at quarter filling. Results are presented for three differ- 
ent final values of J, where the atomic limit (J = 0) 
revival time can be clearly seen to correspond with the 
prediction ij 1 *™ 1 = tt/A. Two effects of finite final hop- 
ping are evident in those plots: first, a clear shift in the 
frequency of the oscillations and, second, a damping of 
the amplitude. In the following, we restrict our analysis 
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to the period and amplitude of the first revival. In the 
homogeneous case, the frequency can be calculated from 
the revival time. In the presence of a confining potential, 
this is in general not possible as the exact revival time 
can change on long time scales. 
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FIG. 1: (Color online) Time evolution of nk=o(t) for three 
final values for the hopping constant, Ja n = 0, Jfl n = 0.3, and 
Jfln = 0.6, in a system with a superlattice potential in ID 
(a) and 2D (b). These calculations were done in a chain with 
L — 400 lattice sites in ID and in a L — 4 x 4 system in 2D, 
both at quarter filling. Time is measured in units of ft/ Jim. 

The quantum corrections due to finite hopping Ai rcv 
versus J for constant A = 1 are presented in a log-log plot 
in Fig. [UJa) for ID and 2D systems. We find from those 
plots that the corrections follow a quadratic behavior for 
values of J < 0.1. This is depicted by a quadratic fit 
f(x) = ax 2 to data points with J < 0.01. Also evident 
from those plots is the very weak dependence of Ai rcv on 
the density both in ID and 2D. This turns out to be very 
convenient later when studying the harmonically trapped 
systems. 

In Fig. HJb) we consider the damping of the oscilla- 
tions, which can be characterized by the amplitude of 
the first revival n T ^ = rik=o(t = t ICV ) subtracted from 
its value in the atomic limit: An™J = n^*°™' rcv — n™I . 
For this quantity we find a quartic behavior, as illustrated 
by the fits in Fig. [2jb) and a much stronger dependence 
on the density. The very fast reduction of the damp- 
ing with decreasing J makes it a less attractive tool for 
experimentally probing small values of J. 

We find numerically the scaling of t rcv with respect 
to the system parameters J and A to have the follow- 
ing functional form t rcv (J, A) = t rcv (J/A)/A whereas for 
the damping n r ^ (J,A) = n' fc c J ( J/A); that is, the for- 
mer depends on the value of A and J/A while the lat- 
ter is only a function of the ratio J/A. In the atomic 
limit, the revival time scales with A: t^° m (A) = tt/A 
and therefore the preceding scaling holds also true for 
Ai rov : At ICV (t,A) = At ICV (t/A)/A. In Sec. HIT] we are 
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FIG. 2: (Color online) Quantum corrections to (a) the revival 
time (Afrcv), and (b) the revival amplitude (AnJ. c Z ), vs J. 
Results are presented for three densities n = 0.125, n = 0.25 
and n = 0.5 in ID and 2D. In (a) a quadratic dependence 
is observed whereas in (b) a quartic dependence is present. 
These behaviors are emphasized by power-law fits for data 
points with J < 0.01. The system sizes for ID are L — 800 
for n = 0.125, and L = 400 for n = 0.25 and n = 0.5, and for 
2D L = 4 x 4 for all densities. Results for n > 0.5 trivially 
follow from the particle-hole symmetry of the model. No data 
are presented for n = 0.125 in (b) because only two particles 
are present in the 4x4 cluster and no damping occurs. 



able to analytically confirm this result for the mean-field 
approximation. On the other hand, in the atomic limit, 
n k=o(A) is independent of A as the system exhibits per- 
fect revivals so Ari™ is only a function of J/A. 



C. Results for hard-core bosons in a trap in ID 

Experimental systems are in general different from the 
ones discussed in Sec. Ill Bl This is because a confining 
potential is required for containing the gas of bosons. The 
confining potential in experiments is to a good approxi- 
mation harmonic, and generates an inhomogeneous den- 
sity profile. Given the results shown in Fig. [UJa), where 
the revival time was shown to depend only weakly on the 
density, one would expect the outcome in the presence 
of a trap not to be strongly dependent on the confining 
potential and the total number of particles. 

Up to small corrections, the preceding turns out to 
be the case for the changes induced in the revival time 
by a finite hopping. However, as shown in Fig. [21a), 
if one quenches J and A keeping constant the trapping 
potential, then a very high damping rate can be seen even 
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in the atomic limit. Hence, measurements at a constant 
curvature of the trap are not the best way to proceed in 
trapped systems. They mix the effects of the trapping 
potential and the finite hopping in the outcome. In fact, 
even the quadratic behavior obvious in the homogeneous 
case (Fig. [2]) becomes obscured in the trap if the confining 
potential is kept the same from the initial state. 
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FIG. 3: (Color online) Time evolution of nk=o(t) in a trap 
after a quench with (a) constant curvature V = 10~ 4 , (b) 
constant characteristic density p — 1.0, and (c) turning off 
the trap (Van = 0), all for the same initial state. In text (b) 
and (c) are referred to as quench type (i) and (ii), respectively. 
These calculations were done for a system with L = 400 and 
N = 100. To keep the characteristic density constant during 
the quench we changed Vim = 10~ 4 — > F nn = for J = 0, 
Fini = KT 4 -> F fl n = 3 x KT 5 for J = 0.3 and F ini = 10~ 4 -> 
Fa,, = 6 x 10~ 5 for J = 0.6. 

In previous work in equilibrium it has been argued that 
the correct way to define the thermodynamic limit for 
a trapped system is by keeping constant the so-called 
characteristic density p = N[V/(dJ)] = , where d is the 
dimensionality of the system (see, e.g., Ref. (I3|)- This is 
equivalent to what is done in homogeneous systems when 
one keeps constant the density N/L. Since in Scc. lIIBl all 
quenches were performed keeping constant N/L, we have 
studied quenches in the trap in which the characteristic 
density is kept constant; that is, one needs to reduce the 
trapping potential by the same amount that the hopping 
parameter is reduced. We denote this scenario quench 
type (i). Another way of reproducing the homogeneous 
results that comes to mind is to remove the trapping 
potential concurrently with the superlattice quench and 
observe oscillations which then take place in a homoge- 
neous potential. This scenario is denoted quench type 
(ii). Within the second approach, one realizes that the 
gas starts expanding after the quench. However, if the 
considered hopping parameters and revival times are suf- 
ficiently small, this will not constitute a problem. Results 
for the dynamics of these cases are shown in Figs. G3R) 
and[3Jc). In contrast to the quench that keeps constant 



the curvature of the trap we now observe that the time 
evolution of nk=o is very similar to the one in homoge- 
neous systems depicted in Fig. [1] 
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FIG. 4: (Color online) Quantum corrections to (a) the revival 
time (At rcv ), and (b) the revival amplitude (An]?Zo) vs J 
for Aa n — 1, both for the quench scenarios (i) keeping the 
characteristic density constant, and (ii) turning off the trap. 
Results are presented for three different initial values of the 
characteristic density p — 1.0, p = 2.0, and p = 3.0, which 
correspond to curvatures Fini = 4.44 x 10~ 5 , Fini = 1.78 x 
10 -4 , and Fni = 4 x 10 -4 , respectively; the system sizes are 
L = 500, L = 400, L = 300, with N = 150 in all cases. Note 
that p — 3.0 already has a Mott-insulating region in the center 
of the trap, while p = 1.0 and p = 2.0 are entirely superfluid. 
The quadratic and quartic fits in (a) and (b) were done for 
points with J < 0.01. 

For both quench scenarios (i) and (ii), we find a 
quadratic behavior for Af rov , which is similar to what was 
shown in Fig. [5] for homogeneous systems. Interestingly 
this behavior is, as depicted in Fig.Uta), practically inde- 
pendent of the quench type and the characteristic density 
of the initial state. We note that for p = 3.0 the initial 
state has an insulating (n = 1) domain in the center of 
the trap, while for the other characteristic densities the 
system is purely superfluid. In ID, the insulator appears 
in the center of the trap when p ~ 2.6 — 2.7 The 
independence of the asymptotic behavior of At rcv on the 
initial state suggests that by measuring the correction to 
the revival time due to the finite value of J in experi- 
ments, it is possible to determine J if one knows A. The 
same can be said for systems without a trap. 

On the other hand, as shown in Fig. |UJb), An™J re- 
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veals a strong dependence on the quench type and the 
initial density. For scenario (i), we obtain a quadratic 
behavior for pure superfluid initial states (p = 1.0 and 
p = 2.0), while for the one having a Mott insulating do- 
main (p = 3.0) a constant damping rate is always present 
for J < 1. For quench type (ii), the damping behaves 
completely differently and shows the quartic behavior ob- 
served in the homogeneous case. With respect to the aim 
to simulate the homogeneous case, this result indicates 
that the fact that the density profile remains unchanged 
during the time evolution under scenario (i) is less im- 
portant than the fact that the potential is homogeneous 
after quench (ii). 

III. MEAN-FIELD APPROACH 



Here ^ = (aj) = Y2=i V^^n-i)"^, = 
a i{n c +i) = 0? an d 2m- denotes summation over all j 
that are nearest neighbors of i. This is a set of L x n c 
equations. The time evolution described by Eq. (|9j pre- 
serves normalization and the total particle number N. 
We solve the system numerically using a fourth-order 
Runge-Kutta method. Self-consistency is guaranteed by 
monitoring the total energy, particle number, and nor- 
malization. 

At this point it is important to stress that this mean- 
field approach is in principle an uncontrolled approxima- 
tion. We gauge its validity against our exact results in 
Sec. IIH 51 Before doing so, we present an instructive ana- 
lytical solution for the equations introduced previously in 
the hard-core limit and for a periodic potential (V = 0). 



Within the mean-field approximation, we can extend 
our analysis to consider the more experimentally relevant 
Bose-Hubbard model: 

H S cb = -J^2(a\aj +H.c.) 

(ij) 

i i 

where [a.;,aj] = Sij and [a^o^] = [aj, at] = 0, as usual 
for bosons. The on-site interaction energy is denoted by 
U. 

The mean-field theory that we employ is based on the 
restriction of the wave function to the Gutzwiller-type 
product state, 



If 



MF/ 



L n c 

i=l 71=0 



(6) 



where n c — > oo for thermodynamic systems, \n)i denotes 
a single-site Fock state for lattice site i and the complex 
coefficients a.i n allow for a time dependence. For all nu- 
merical calculations, a finite cutoff n c is taken. 

The mean-field ground state in equilibrium is found by 
minimization of the energy expectation value, 



(^mfI^scb -fJ.N\V 



MF 



(7) 



where p is the chemical potential and N counts the total 
number of particles. Hence, from here on we work on 
the grand-canonical ensemble. To find the time-evolution 
of the mean-field approximated system, we employ the 
time-dependent variational principle [18| that minimizes 
the expression 

(*Mf|*$ -ffscB +pN\^ MF ) , (8) 
and yields the following set of differential equations: 

iain = ^ (Vn + laj(„ + i)$* + \^aj( n -i)$j) 



+ OLin n 



u 



(n - 1) + Vrf - p 



(9) 



A. Analytical mean-field solution for hard-core 
bosons in a periodic potential 

In the case of hard-core bosons, it is possible to re- 
duce the number of equations considerably and employ 
a parametrization for the cti n in Eq. §§§ that preserves 
normalization and deals with real variables - this is due 
to the equivalence of hard-core bosons to s = 1/2 spins, 
which leads to the following ansatz for the Gutzwillcr 
wave function I19i: 
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! = 1 



sin — 
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cos — e 



^4 |0>. (10) 



If there is no trap in the system (V = 0) it is possible 
to use translational invariance to simplify the equations 
(0), in which in the hard-core limit we again employ the 
superlattice quench introduced before. This leads to a 
formal replacement of U{n— 1)/2 by A in ([9]). As all sites 
with the same potential must have the same properties 
and the system's wavefunction is a product of single-site 
states, the L-site system can be reduced to an effective 
two-site problem (for two on-site potentials p\l2 = M^^) 
independent of the dimension. Then the ansatz Eq. (|10[) 
yields the following form of Eq. © : 

B x = -2d Jsin6» 2 sin0, (11a) 
2 = 2d Jsin^isin^, (lib) 
4> = 2A — 2 d J(sm02 cot#i — sin^i cot #2) cos</>, (11c) 

where <f> = (f>\ — <f>2- Here it can be seen that dimension- 
ality enters the equations only by a simple rescaling of 
the hopping parameter: J — > d J. 

Also, the argument of translational invariance allows 
one to find a simple expression for n^o in the two-site 
system: 



rifc=o = n + j sin B\ sin 82 cos < 



(12) 



Figure O depicts the time evolution of Uk=o for the same 
systems for which the exact solution was presented in 
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FIG. 5: (Color online) Plot of the time evolution of hard-core 
bosons in a periodic potential in the mean-field approximation 
for ID (a) and 2D (b) after a superlattice quench = — > 
.Afin = 1 for a system at quarter filling. This is to be compared 
with the exact results in Fig.[TJ where exactly the same system 
parameters were used. 

Fig. [5J One can clearly see that the mean-field and exact 
results show a similar shift of the frequency. However, 
the mean-field solutions exhibit no damping, and hence 
they are qualitatively incorrect for that quantity. 

At this point it is instructive to extract the result for 
the atomic limit: for J = 0, Eq. ([TTj) has the trivial 
solution (f)(t) = 2At + </>(0) and #1/2 constant. Insertion 
of this result in Eq. (|12[) immediately reveals the revival 
time tf c ° m = n/A. 

Obtaining the solution of the system of Eqs. (fTTj) for fi- 
nite J is possible by treating them like a classical system. 
Identification of Hamilton functions and the observation 
of its trajectories leads to an analytical expression for the 
period of rik=o: 

ticv = / du/(u), where (13) 

/(«) - {d 2 J 2 (l - u 2 )[l - (2 7 - uf] - (H' - Auf}^ 

where 7 = 2n — 1 and H' Q = — 4rt(l — n) dJ—jA. A closed 
expression for the preceding integral exists. However, it is 
cumbersome and does not provide any apparent informa- 
tion on the functional form of t Tev as it depends on the el- 
liptic integral of the first kind. The integral limits u\ and 
U2 are the solutions of 1/ f(u) = that lie within [— 1, 1]. 
This requires solving the root of a polynomial of fourth 
order, which can also be done analytically. The lower 
boundary u\ is given by a simple expression: u\ = 7. In 
the case of half filling, also the upper boundary is given 

by a simple expression: u 2 = ^jj + 8d 2 J 2 /A 2 — 1J . 

In Fig. [6l we plot the analytic solution for different 
dimensions at quarter filling. As mentioned before, di- 
mensionality in the mean-field picture is captured by a 
simple rcscaling of J — > d J. For comparison, we de- 




FIG. 6: (Color online) Plot of the revival time At rev = t^° m - 
irev for a system of hard-core bosons in a periodic potential 
and a superlattice quench A- m i — — > Afi n — 1. t rov is given 
by the analytic solution Eq. (|13|l . Depicted are results for ID, 
2D, and 3D. The data points are numerical solutions of the 
set of Eqs. © . 

pict numerical solutions of Eq. (|9|) as points in the plot. 
The latter are required for studying the inhomogeneous 
trapped hard-core boson case and soft-core bosons, for 
which no analytic solutions are available. 

The analytic expression (|13|) allows us to confirm the 
numerical finding for the scaling relations of t rev and 
with respect to the parameters A and J. The proposed 
scaling for the revival time t TCV (J, A) = t Iev (J/A)/A does 
obviously hold for the integrand f(u) in Eq. (fT3|) . The 
fact that for the integration limits one has u 1 / 2 (J, A) = 
u±/2(J/A), proves the proposition for the mean-field case. 

Furthermore, the evaluation of Eq. (fT2"j) leads to the 
expression 

n k=0 =n- 4Zj(H' q + A cos 0i) , (14) 

where one can see, (i) that no damping of nk=o occurs 
within the mean-field approximation in the hard-core 
limit and (ii) that n k =o only depends on J/ A as found 
for the numerical solution. 

We note also that for the case of half filling and d J = 1, 
Eq. ([TB")) yields t ICV — > 00, reflecting the fact that the 
mean-field equations of motion do not predict any oscil- 
lations in this case - in contrast to the exact solution. 
With the absence of damping and the last observation, 
there are already two deficiencies of the mean-field ap- 
proximation that we need to keep in mind for the analy- 
sis that follows - this makes the comparison to the exact 
solution an essential duty to ensure one has an idea of the 
limits of the validity of the mean-field results presented 
in Sec. El 



B. Exact vs mean-field results 

In equilibrium, a detailed comparison between the pre- 
dictions of the mean-field theory introduced before and 
exact quantum Monte Carlo simulations for the ground 
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state of hard-core bosons in the presence of a superlattice 
potential was presented in Refs. [13, El- The Gutzwiller 
approach was found to correctly predict the two phases 
present in the ground state of this model, namely, a su- 
pcrfiuid phase for all fillings but n = 0, 1/2, and 1 and for 
n = 1/2 below a critical value of A/ J and a Mott insula- 
tor (a charge density wave) for n = 1/2 above a critical 
value of A/ J. However, Gutzwiller mean- field theory was 
shown to provide a poor estimate of the critical value of 
A/ J for the superfluid-Mott-insulator transition. It over- 
estimated it by around 100% in 2D and a 50% in 3D. 

As shown in the previous section, after the quench in 
the periodic system, the mean-field solution does not ex- 
hibit any damping for the amplitude of the oscillation 
whereas in the exact solution there obviously is damp- 
ing. For this reason, we do not study the damping any 
further. In the remainder of the article we therefore focus 
on the predictions of mean- field theory for the corrections 
to the revival time. 

In order to be more quantitative, wc define the relative 
deviation of the mean-field approximation from the exact 
solution by 



= At- (J) - At?J v (J) 



A+cx 
'-"'rev 



(-0 



(15) 



where Ai°* v 



(J) and At™ { v (J) are the corrections to the 
revival time due to a finite value of J for the exact and 
mean-field solutions, respectively. 
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FIG. 7: (Color online) Comparison of the mean-field approxi- 
mation and the exact solution for ID (a, c), and 2D (b, d). In 
(a) and (b), we plot the error introduced in (|15|) on a linear- 
log scale for the systems presented in Fig. [2] In (c) and (d), 
we again show the results presented in Fig.[2]but this time on 
a linear scale together with the mean-field results, where the 
latter are drawn as lines and were calculated via the analytical 
solution (TT5)| . 

In Figs.JTJa) and[7|b), we plot e( J) in a linear- log scale 
vs J. That figure shows an almost constant error over 



two decades (10 -3 < J < 10" 1 ). For J < 10~ 3 rounding- 
off errors set in as one starts dealing with numbers ~ 
10~ 10 ; that is, the deviations seen in that region are not 
to be considered any further. The absolute value in the 
constant part of the deviation is around 16% for ID (a) 
and 11% for 2D (b) for n = 0.125 and n = 0.25. At half 
filling, interestingly, the deviation is yet much smaller. In 
all cases, it is obvious that the mean-field approximation 
describes the 2D system better than the ID case, even 
though the system size in ID is one order of magnitude 
larger than in 2D and mean-field is expected to be more 
accurate as the system size is increased. 

In Figs.[7](c) and[7{d), we present the same results as in 
Fig. [3 and compare them to the mean-field predictions, 
but this time in a linear scale. This scale emphasizes 
the differences between the mean-field results and the 
exact ones for values of J close to A = 1. Once again, it 
is obvious that the mean-field approximation works bet- 
ter in 2D than in ID, and that it becomes a very good 
approximation of the exact results for n = 0.125 and 
n = 0.25. Wc further note the already-mentioned case of 
d J = 1 , which does not yield any revival in the mean-field 
approximation and therefore Ai rcv — > oo as the figures 
suggests. Interestingly, there is a corresponding anomaly 
in the exact solution in 2D: there, in a neighborhood of 
J = 0.5, the exact solution does not produce a symmet- 
ric peak after the first revival oscillation, which makes 
it meaningless to determine a value for the revival time 
- therefore these data points are missing. Finally, Figs. 
E[a) and (b) show that at half filling mean-field theory 
provides the most accurate results for the quadratic re- 
gion of small values of J, whereas it provides the worse 
results for large values of J, as seen in (c) and (d). 




10"* 10" 3 10" 2 10" 1 10° 10" 4 10" 3 10" 2 10" 1 10° 
J J 

FIG. 8: (Color online) Plots of the error ([I5j| for the ID 
trapped case presented in Fig.[4]for: (a) quench type (i) keep- 
ing the characteristic density constant and (b) quench type 
(ii) turning off the trap, both presented on a linear-log scale. 

Results for e(J) in the ID trapped case are presented 
in Fig. [5] for the two types of quenches introduced in the 
previous section: (i) the scenario of a constant character- 
istic density and (ii) switching off the trapping potential. 
The behavior is qualitatively the same as just discussed 
for the homogeneous case of Fig. EJa). Quantitatively, 
we find an error ~ 5%, which is in between the values for 
low and half filling in Fig. [TJa), and is similar for both 
quench types. Such an intermediate value is expected 
because the trap causes a density profile with different 
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densities in different regions of the trap. 

For the description of experiments, which are of more 
interest in 3D trapped geometries and very large system 
sizes, one can expect much smaller errors than the ones in 
Fig.[5J Therefore, we conclude that the shift of the revival 
time due to a finite hopping is correctly captured not 
only qualitatively but also quantitatively by the mean- 
field approximation described here. This is an interesting 
finding considering that, in contrast, the description of 
the evolution of the amplitude in terms of the mean-field 
approximation is incorrect even at the qualitative level. 



IV. RESULTS FOR SOFT-CORE BOSONS 

In Refs. [H, [2(| it has already been shown that ultra- 
cold bosonic gases in optical lattices can be well described 
by the bosonic Hubbard model (0. In light of the re- 
cent results in Refs. [fl Q mentioned in the introduc- 
tion, we note that the interaction constant U and the 
hopping amplitude J in the Hamiltonian ([5]) are the ef- 
fective two-body interaction and one-body hopping, re- 
spectively, whose origin is not under discussion here. 
They may be obtained after multi-orbital rcnormaliza- 
tion effects are taken into account. Multi-orbital effects 
may also generate effective higher-body interactions that 
translate into additional frequencies during the collapse 
and revival of the matter-wave interference but are not 
considered here. These effects can be reduced by prop- 
erly engineering the initial state. We focus on the effect 
that a finite effective hopping J has on the fist revival of 
the matter wave. 

Collapse and revival oscillations like the ones observed 
experimentally in Refs. [l], H[ nave been reproduced in 
ID bosonic 21 1 and fermionic [22j systems by means of 
numerically exact time-dependent renormalization group 
techniques. Here we use mean-field theory to study 3D 
bosonic systems theoretically. Following the results of the 
previous section, we expect the mean-field predictions for 
the tunneling- induced correction of the revival time of the 
value in the atomic limit Ai rev to be close to the exact 
results. 

We proceed in the way introduced in the preceding 
section. For soft-core bosons, we solve the equations 
© for a cutoff of n c = 7 numerically. This cutoff en- 
sures convergence with respect to the energy expectation 
value and the momentum distribution, and was success- 
fully employed in Ref. [23| to study the superfiuid-Mott- 
insulator transition. Initially, the system is prepared in 
the Gutzwillcr mean-field ground state of the trapped 
system at an intermediate interaction U- m i = 6 Ji n ; ( J; n ; 
1 sets our energy scale), which ensures the validity of a 
one-band model in experiments while the system is still 
far from the transition to the Mott insulator. Within the 
mean-field approximation ([// J) cr it = 34.8 (for six near- 
est neighbors) [H, 0, IH- At time 1 = °- tnc on-site 
interaction is doubled to [/fi n = 12 and we investigate 
the collapse and revivals for several values of Jfi n < 1 



(where, as in the section about hard-core bosons, the no- 
tation [/fin = U and Jfi n = J is used in all unambiguous 
cases) . 




FIG. 9: (Color online) Dynamics of the zero momentum peak 
nk-o vs time t for a system with homogeneous potential and 
density n = 1.5. Time is measured in units of ft/Jini. These 
results are derived from the mean-field solution of the time 
evolution ([9]) for the Bose-Hubbard model after an interaction 
quench from U- m i = 6 to U^ n = 12. Several final values of the 
hopping constant J are depicted. 

Figure[9]depicts the collapse and revival oscillations for 
the homogeneous case. These predictions, within mean- 
field theory correspond to the solution of a single site 
problem. They clearly exhibit a different functional form 
when compared to the ones for hard-core bosons in a 
superlattice potential Fig. [2j Also, not only a shift in 
the revival time is present but additionally a damping 
of the amplitude can be observed. We do not consider 
this damping effect, as we cannot make any statement 
about the validity of the mean-field approximation for 
this quantity (see Sec. IHIj) . We note again that the re- 
vival time for the interaction quench in the atomic limit 
is given by £ rcv = 2-k/U. 

In Fig. QJ2 we show At rGV for several densities in the 
homogeneous case. We observe a linear relation empha- 
sized by the fits for data points with J < 0.01 in the 
figure. Since soft-core bosons are not subject to particle- 
hole symmetry, the behavior with increasing density is 
different from the one observed for hard-core bosons in 
Fig. [2j The smallest deviation from the atomic limit is no 
longer reached at half filling; instead, for soft-core bosons 
we find that it appears for a density n ~ 0.75. Except 
for densities around this value, Ai rev is not strongly de- 
pendent on the density. With respect to the dependence 
of Ai rev on both parameters J and U, we note that the 
same scaling as in the hard-core case (with A and U in- 
terchanged) holds true: t iev (J,U) = t ICV (J/U)/U. 

The calculations for the trapped case are more de- 
manding computationally. This is because translational 
invariance is broken and one has to deal with all the lat- 
tice sites in the system. The results reported here are 
obtained for a lattice with L = 30 x 30 x 30 = 27000 sites 
with N = 1000 up to N = 11000. As before, finite size 
effects for our observables of interest are extremely small. 
As a matter of fact, we found that it would be difficult 
to distinguish the results reported here from those of a 
L = 10 3 = 1000 system. Again, results are presented for 
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FIG. 10: (Color online) Plot of At ICV vs J for soft-core bosons 
in a homogeneous system and different densities after an in- 
teraction quench from f/jni = 6 to [/ nn = 12. This is done for 
a single-site system, as the homogeneous result is independent 
of the system size in the mean-field approximation. Note that 
the results for different densities lie very close to each other. 
The linear fits are done for data points with J < 0.01. 

the two quench scenarios analyzed in detail in Sec. Ill CI 
for hard-core bosons. 
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FIG. 11: (Color online) Plot of the time evolution of the 
momentum peak rik=o for ~ 11000 soft-core bosons in a lattice 
with L — 30 x 30 x 30 sites for: (a) quench scenario (i), and 
(b) quench scenario (ii) (see text). The initial state has a 
characteristic density of p = 39.8 and an on-site energy LAni = 
6. At time t = the interaction is quenched to Uf\ n = 12 and 
the curvature of the trap is modified according to (i) and (ii). 

Results for the time evolution of rik=o m the harmonic 
trap are shown in Fig. [TTJ for an initial state with char- 
acteristic density p = 39.8 and the two different quench 
types: (i) keeping constant the characteristic density and 
(ii) turning of the trap. Here, we observe an effect that is 
qualitatively different from the the one seen in the hard- 
core limit (Fig. [3]) and the homogeneous soft-core case, 
namely, the revival time in the case with finite hopping 



exceeds the atomic value [Fig. [Til a)]. This effect is only 
observed in the quench scenario (i) for high characteristic 
densities. For quench type (ii) the effect is not present 
for any density [Fig. [TTTb)]. 
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FIG. 12: (Color online) Quantum corrections to the revival 
time (At ICV ) vs J for the scenarios: (a) keeping the charac- 
teristic density constant, and (b) turning of the trap; and 
an interaction quench [/mi = 6 — ¥ Ufi n — 12. Results are 
presented for five different initial values of the characteristic 
density - p = 0.6, p = 3.1, p = 10.6, p = 23.8, and p = 39.8 - 
that correspond to the following values of the density in the 
middle of the trap: n ce nter ~ 0.5, reenter ~ 1.0, ^center ^ l.O, 
"■center ~ 2.0 and n C enter ~ 2.5. Note that in (a) At rcv has 
a negative sign for p = 23.8 and p = 39.8; therefore, the 
absolute value is depicted. 

In Fig. [T21 At rov is plotted for several values of the 
characteristic density of the initial state. In Fig. IT^T a). 
results are shown for quench type (i). Note that we ac- 
count for the fact that the revival time exceeds the atomic 
limit for high densities by plotting the absolute value of 
At rcv ; in particular, p = 23.8 and p = 39.8 yield a neg- 
ative value of At r ev The strong dependence of the re- 
sults on the characteristic density of the initial system, or 
equivalently, the initial density in the center of the trap, 
makes this scenario unsuitable for experimentally prob- 
ing the values of J after the quench. The experimental 
uncertainty of the filling in the center of the trap would 
lead to a large uncertainty in determining J. 

Scenario (ii) seems to be a good candidate for the latter 
goal. As depicted in Fig. [T2Tb), for characteristic densi- 
ties that are not too small (p > 10) i.e. for densities in 
the center of the trap that are ^center ^ 1-5, we observe 
a weak dependence of the revival time on the character- 
istic density of the initial state. This is usually fulfilled 
in experiments like [l[. We also note that the relative 
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deviation At' rcv = At Iev /tfl° m is only one order of mag- 
nitude smaller than the normalized hopping parameter 
J/U, due to the linearity of the relation and a prefactor 
~ 0.1. This shows that the described effect is not small 
and one should be able to measure it in experiments. 

Finally, we should mention that we also performed cal- 
culations for different values of the interaction constant U 
before and after the quench. They all exhibited a similar 
qualitative behavior as depicted in Fig. [12j We there- 
fore stress that our results do not depend on a particular 
value of U but represent a general behavior that can be 
reproduced with experimentally relevant parameters. By 
comparing experimental results and calculations within 
mean-field theory, plus using the expected experimental 
values for U, one could then use experimentally measured 
values of Ai rev to determine the final value of J. 

V. CONCLUSION 

We have presented a detailed study of the dependence 
of collapse and revival properties of the matter-wave in- 
terference in lattice boson models on a finite tunneling 
amplitude after an interaction quench. 

We first studied quenches of hard-core bosons on a su- 
perlattice potential. For those systems, we presented 
exact numerical results in ID and 2D, and compared 
them with the approximated mean-field solution. Both 
approaches exhibited the same functional form of the 
correction to the revival time produced by finite hop- 



ping parameters after the quench, with a leading behav- 
ior ~ t 2 /A 3 . The mean-field results were also shown to 
have, as expected, smaller errors in 2D as compared to 
ID. Since the largest errors for homogeneous 2D systems 
were ~ 10% and for trapped ID systems ~ 5% (in con- 
trast to ~ 17% for the ID homogeneous case), we expect 
that in 3D trapped systems, mean-field theory should 
provide relatively accurate results for the corrections to 
the revival time. 

We then studied interaction quenches in the Bose- 
Hubbard model in 3D. In this case, our analysis was 
solely based on the Gutzwiller mean-field theory. We 
showed that for soft-core bosons the corrections to the 
revival time in the atomic limit, due to finite values of 
J after the quench, are ~ J/U 2 . This is an effect that 
could be measured experimentally. Given the weak de- 
pendence of the correction on the initial density profile, 
provided the density in the center of the trap is greater 
than n = 1.5, we have proposed that the corrections to 
the revival time measured experimentally could be used 
to determine the actual value of J after the quench. The 
only input one would need is the experimental value of U 
and the mean-field predictions from calculations similar 
to the ones presented here. 
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